import numpy as np
from scipy.stats import norm

# 清除变量（Python 中无需显式清除，但可初始化）
# 参数赋值
T = 1        # 到期时间（年）
r = 0.04     # 无风险利率
sigma = 0.3  # 波动率
X = 5        # 执行价格（行权价）
S = 6        # 当前股票价格

# 计算 Black-Scholes 显式解（解析解）
# d1 和 d2 是 Black-Scholes 公式中的关键参数
d1 = (np.log(S / X) + (r + 0.5 * sigma**2) * T) / (sigma * np.sqrt(T))
d2 = d1 - sigma * np.sqrt(T)

# 欧式看涨期权的解析价格
V_ex = S * norm.cdf(d1) - X * np.exp(-r * T) * norm.cdf(d2)
print(f"解析解期权价格: {V_ex:.6f}")

# 模拟数量（样本数）
N = 1000000

# 生成标准正态分布的随机数（相当于 MATLAB 的 randn）
Z = np.random.standard_normal(N)

# 根据几何布朗运动模型生成到期日股票价格样本
# S_T = S * exp((r - 0.5*sigma^2)*T + sigma*sqrt(T)*Z)
S_T = S * np.exp((r - 0.5 * sigma**2) * T + sigma * np.sqrt(T) * Z)

# 计算每条路径下的期权收益（payoff），贴现后取期望
# payoff = max(S_T - X, 0)，然后乘以 e^(-r*T)
V_mc = np.exp(-r * T) * np.maximum(S_T - X, 0)

# 计算样本均值（即蒙特卡洛估计的期权价格）
V_mc_mean = np.mean(V_mc)
print(f"蒙特卡洛估计期权价格: {V_mc_mean:.6f}")

# 计算绝对误差和相对误差（误差点）
err = abs(V_mc_mean - V_ex) / V_ex
print(f"相对误差: {err:.6f}")

